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ABSTRACT 

Context. With the plethora of detailed results from heliospheric missions such as Ulysses and SOHO and in advent of the first mission 
dedicated to in situ studies of neutral heliospheric atoms IBEX we have entered the era of precision heliospheric study. Interpretation 
of these data require precision modeling, with second-order effects quantitatively taken into account. 

Aims. We study the influence of the non-flat shape of the solar Lyman-o- line on the distribution of neutral interstellar hydrogen in the 
inner heliosphere and assess importance of this effect for interpretation of heliospheric in situ measurements. 

Methods. Based on available data, construct a model of evolution of the solar Lyman-ct line profile with solar activity. Modify an 
existing test-particle code calculating distribution of neutral interstellar hydrogen in the inner heliosphere to take into account the 
dependence of radiation pressure on radial velocity. 

Results. Discrepancies between the classical and Doppler models appear at ~ 5 AU and increase towards the Sun from a few percent 
to a factor of 1.5 at 1 AU. The classical model overestimates density everywhere except a ~ 60° cone around the downwind direction, 
where a density deficit appears. The magnitude of discrepancies depends appreciably on the phase of solar cycle, but only weakly 
on the parameters of the gas at the termination shock. For in situ measurements of neutral atoms performed at ~ 1 AU, as those 
planned for IBEX, the Doppler correction will need to be taken into account, because the modifications include both the magnitude 
and direction of the local flux by a few km/s and degree, which, when unaccounted for, would bring an error of a few degrees and a 
few km/s in determination of the bulk velocity vector at the termination shock. 

Conclusions. The Doppler correction is appreciable for in situ observations of neutral H populations and their derivatives performed 
a few AU from the Sun. 

Key words. 



1. Introduction 



After discovery at the end of 1960- ties of a diffuse interplanetary 
k* - ■ Lyma n-o' glow (Thomas & Krassal fl97 it Bertaux & Blamont 
1 1197 lb . predicted bv lFahrf d 1968b and lBlum & Farm d 1970b as due 
, to scattering of solar Lyman-ff radiation on the neutral interstel- 
' lar gas flowing through the Solar System, development of mod- 
els of the distribution of this gas in the heliosphere begun. 

At the very early phase, the influence on neutral heliospheric 
gas of processes going on at the heliospheric interface was ne- 
glected. It was assumed that both the ionization due to solar out- 
put and the radiation pressure acting on the H atoms are station- 
ary and spherically symmetric around the Sun and that they fall 
off proportionally to inverse square of heliocentric distance. It 
was further assumed that the inflowing neutral gas is monoen- 
ergetic, i.e., that before the encounter with the Sun the atoms 
move with identical, parallel-oriented velocities, identical with 
the macroscopic bulk velocity of the gas. 

These assumptions formed basis of the first model of density 
distribution of neutral int erstellar gas around the Su n: the purely 
analytical "cold" model dFahj|l968t lAxfordlll972l) . It features 
axial symmetry about the axis of inflow of the gas and shows ei- 
ther a singularity at the downwind axis, when radiation pressure 
is too weak to compensate solar gravity, or an empty "avoidance 
zone" in the downwind region with a paraboloidal boundary sur- 
face, when radiation pressure overcompensates solar gravity. 
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Lifting of the monoenergetic assumption (i.e. allowing for a 
finite temperature of interstellar gas, high enough to yield ther- 
mal velocity comparable to the bulk velocity of the gas), adop- 
tion of the distribution function of the gas far away from the Sun 
("in infinity") in the form of a Maxwellian shifted in the velocity 
space by the bulk velocity vector, and assumption that the gas is 
collisionless on the distance scale comparable to the size of the 
heliosphere allowed to use the Boltzmann equation to describe 
the problem of interaction of neutral interstellar gas with solar 
environment. Its solution brought t he "hot model" of the gas 
distribution in the inner heliosphere dThomas|[l97^ :lFahiiri978. 
ll979tlWu & Judgdl979tlLallement et al.ll985al) . This model re- 
quired numerical integration of the distribution function, but the 
function itself was still analytical. The "hot model" in its vari- 
ous implementations became the canonical model of neutral in- 
terstellar gas distribution in the inner heliosphere. 

Further development of modeling of the interaction and dis- 
tribution of neutral interstellar hydrogen near the Sun focused on 
two main topics. On one hand, a lot of effort was put to under- 
stand and simulate processes going on at the boundary between 
the expanding solar wind and the incoming partially ionized in- 
terstellar gas, in the region referred to as the heliospheric inter- 
face - and this issue will not be addressed i n the pre s ent pap er: 
the reader is referred to recent reviews by lBaranovl d2006bllah : 
Ilzmodenovid2006h ; llzmodenov & Bara nov (2006). On the other 
hand, development of the hot model continued, aimed at a more 
realistic description of distribution of neutral interstellar hydro- 
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gen in the inner heliosphere, suitable for quantitative, and not 
only qualitative interpretation of heliospheric measurements. 

In the first shot, the assumption that the i onization rate is 
spheri cally symmteric was eliminated, when Lallemen t et al.l 
(1985b) described the latitudinal modulation of the charge ex- 
change rate with a one-parameter formula 1 - A sin 2 <f>, imple- 
menting it in the CNRS model of the heliospheric gas distri- 
bution. This allowed to vary the equator-to-pole contrast of the 
ionization rate, but required to keep fixed the width and range of 
the enhanced ionization band. 

A different extension of the hot model was proposed by 
IRucinski & Fahl d!989lfl99lh . who realized that the rate of ion- 
ization by electron impact is not proportional to inverse square 
of solar distance. Electr on ionization is of particular importance 
for interstellar he lium (McMulli n et al.l l2004t ILallemen t et alj 
l2004t Iw"ittell2004l) ; for hydrogen it is noticeable inside a few AU 
from the Sun, where the density of hydrogen gas is already very 
much reduced by earlier ionization and solar radiation pressure 
overcompensating solar gravity. Therefore this aspect of physics 
of neutral interstellar hydrogen was neglected for quite a while 
since t hen, but is reintroduc ed in the most recent versions of the 
model dBzowski et al.ll2008l) . 

In the next round of development of the density model, 
implemented in the Warsaw test-particle code, further on 
referred to as the Warsaw model, gone was the assump- 
tion of invariability of radiation pressure and ionizatio n 
rate dRucinski & Bzowskil 1995at IBzowski & Rucinskilll995allbl 
IBzowski et al.l 1 19971) .Indeed, both th e solar EUV flux 
dFlovd et al.l l2002) and the solar wind flux (King & Papit ashyiiil 
120051) vary considerably during the solar cycle. The result is a 
solar cycle variation of the solar radiation pressure, of the EUV 
ionization rate, and of the rate of charge exchange between neu- 
tral H atoms and solar wind protons. This in turn results in ap- 
preciable variations of density and bulk velocity of neutral inter- 
stellar hydrogen within a dozen of AU from the Sun. 

The Warsaw time dependent model was fully numerical be- 
cause not only the integration of the distribution function had to 
be performed numerically, as in the classical "hot model", but 
also the calculation of the integrand, i.e. of the distribution func- 
tion itself. At this phase of heliospheric research, in lack of suf- 
ficiently long time series of measurements of solar wind speed 
and density, and of the solar EUV output, substitute idealized 
models of evolution of these parameters had to be used. 

These aspects of development of heliospheric gas mod- 
els were discussed in greater detail in a review by 
IRucinski & Bzowskil d!996l) . 

In the next move, effects of interaction of the solar wind and 
interstellar gas in the heliospheric interface were taken into ac- 
count. Owing to the charge exchange between the atoms of in- 
terstellar gas and the heated and compressed plasma in front of 
the heliopause, the original population of neutral atoms is some- 
what cooled and accelerated, and a new population of atoms 
appears. These new neutral atoms originally inherit the proper- 
ties of the parent plasma population, but further decouple from 
this plasma and flow through the heliopause and inside the ter- 
mination shock of the solar winds. Both populations interact 
further with the ionized components, exchanging charge with 
the protons from local plasma. Hence, the processes in the he- 
liospheric interface crea te a few distinct , collis i onless popula- 
tions o f neutral atoms dOsterbart & Fahrl dl992f) ; iBaranov et alj 
(199lbjB aranov & Malama ( 1993]), as dis cussed in detail by 
Ilzmodenovld2000h and lMalama et all d2006ln . 

From the view point of modeling of the distribution of neu- 
tral interstellar hydrogen in the inner heliosphere, the most im- 



portant aspect of the processes going on in the heliospheric inter- 
face is the modification of distribution function at TS. Instead of 
the shifted Maxwell ian with parameters homogeneous in space, 
IScherer et al.l (1999) adopted a sum of two Maxwellians, with 
non-isotropic temperatures, shifted by appropriate bulk velocity 
vectors. The two components of the new functions had parame- 
ters being functions of the offset angle 6 from the upwind direc- 
tion and corresponded to the two thermal populations (primary 
and secondary), as predicted by the Moscow Monte Carlo simu- 
lation of heliospheric interface. This new version of the Warsaw 
model was time dependent, but still axially symmetric. At this 
time, enough measurements had been published to attempt to 
introduce observations-based models of radiation pressure and 
ionization rate to the simulations. 

Axial symmetry was removed in the next round of the 
Warsaw model develop ment, when anisotropy of th e ioniza- 
tion rate was introduced. IBzowski et a D d200lL |2002) allowed 
the ionization rate to change as a continuous function of heli- 
ographic latitude, with the latitudinal profile of the ionization 
rate continuously changing with the phase of solar cycle. Hence 
the ionization field in the model became 2D (keeping the axial 
symmetry about the solar rotation axis) and time-dependent, and 
since the gas inflow axis and solar rotation axis are inclined at an 
angle to each other, the model of neutral hydrogen distribution 
in the inner heliosphere became 3D and time dependent. 

The most recent extension of the Warsaw model is presented 
in this paper. We improve on the modeling of radiation pressure 
acting on individual H atoms, which now is not only a function 
of time, but also of radial velocity of the atom. In the previous 
versions of the model it was assumed that the profile of the solar 
Lyman-ff line, responsible for the radiation pressure, is flat. In 
reality it not only is non-flat, but shows co nsiderable variations 
depending on the p hase of solar activity dLemaire et alj f2002; 
iLemaire et al.ll2005l) . 

We take this into account and in Section 2 we develop an 
observation-based model of evolution of the line profile as a 
function of the line- and disk-integrated flux. We use this model 
in a newly-developed code which simulates density and higher 
moments of distribution function of neutral interstellar hydrogen 
in the inner heliosphere. We discuss these calculations in Section 
3. With the use of the newly-developed code we assess modifi- 
cations of the local gas density and of its local flux with respect 
to the results of models neglecting the dependence of radiation 
pressure on radial velocity. We also show possible implications 
for interpretation of heliospheric in situ measurements of neutral 
atoms and pickup ions. This discussion is provided in Section 4. 
Section 5 offers a summary of the results. 



2. The Doppler model of neutral hydrogen 
distribution in the inner heliosphere 

The approach exercised in the present model is a modifica - 
tion of the approach presented by Ruc inski & Bzowskil l[l995a). 
Density and higher moments of distribution function of neutral 
interstellar hydrogen at a location R in the inner heliosphere at a 
time r are calculated by numerical integration of the local distri- 
bution function constructed as a product of the distribution func- 
tion in the source region w src and of the probability w loa of sur- 
vival of the atoms traveling from the source region to the local 
point R against ionization: 

/ (R, V, t) = w src (r SK (R, V, t) , v src (R, V, t)) w ion (R, V,r) (1) 
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Fig. 1. xa Selected profiles of the solar Lyman-o' line observed 
by SUMER/SOHO (dots) for solar minimum (lower line, Feb. 
17, 1997) and solar maximum (upper line, May 20, 2000), com- 
pared with the results of the fitted model, specified in Eq. (O 
(lines). The gray part of the model lines corresponds to ±30 km/s 
around and illustrates the range on the profile which is the most 
relevant for the thermal populations of neutral H atoms in the 
heliosphere. The horizontal lines mark the yU values used in the 
comparison simulations with the use of the classical hot model. 
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Fig. 2. xb Comparison of models of the solar Lyman-o' line pro- 
file: solid line: the present model, dashed line: Chabrillat & 
Kockart, dashed-dot line: Fahr et al.; dots: Scherer. 



where r stc (R,V,T) is the start position in the source region of 
the atom that at the local point R at time t has velocity V, and 
v src (R, V, t) is relevant start velocity in the source region. 

What concerns the code, the distribution function in the 
source region w slc a priori can be any reasonable function. In 
the following we adopted it either as a Maxwellian shifted by 
the bulk velocity vector, homogeneous in space and independent 
of time: 



WW (r src , v src ) = n s 



2nkT*, 



3/2 



exp 



m (V B - Vsrc)" 

2.k Terr 



(2) 



or, in the two-populations case, as a sum of two such func- 
tions with different parameters n slc (density in the source re- 
gion), Vb (bulk velocity at the source region), and T src (tem- 
perature at the source region); m is the atom mass, and k the 
Boltzmann constant. Basically, the parameters of these popula- 
tions are functions of th e location in the heliospheric interface 
dlzmodenov et al. 2001). Since our goal in this paper is investi- 
gating the effect of the v, -dependence of radiation pressure on 
the gas distribution, we left this problem out of the analysis. 
Some insight into modifications of the gas density profile due 



to the angular gradient of the parame ters of the H populati ons at 
the termination shock is provided bv lBzowski et al.l (|2008). 

The challenging issue in the calculation of the local distribu- 
tion function of the gas is finding of the link between the local 
velocity of the test atom V at time r and location R and its po- 
sition r src and velocity v src in the source region in the scenario 
where the radiation pressure is a function of time and radial ve- 
locity. Indeed, the force acting on a H atom in the heliosphere is 
composed of the attracting solar gravity and the repelling radia- 
tion pressure, which is a function of the spectral flux 1^ relevant 
for a given radial velocity v r . The spectral flux 1^ is also a func- 
tion of the line-integrated flux 7 tot (f), which varies on time scales 
from days to dozens of years. Thus the spectral flux responsible 
for instantaneous radiation pressure acting on a H atom moving 
with radial velocity v r is a function of v r , / to t (0 and, implicitly, 
of t:h = h(y r ,I tot (t)). 

Since both 7 tot and solar gravity fall off with heliocentric dis- 
tance as l/r 2 , the radiation pressure can be expressed by the fac- 
tor fi of compensation of solar gravity. In the modeling, this fac- 
tor is a function of radial velocity and of the line-integrated flux, 
which in turn is a function of time. Hence finally fi = fi(v r , t). 

The link between the local position and velocity vectors R, V 
and the corresponding position and velocity vectors at the source 
region r src , v src is determined by numerical solving of the equa- 
tion of motion specified in Eq.Q, performed backwards in time. 
Along with the calculation of the trajectory of the atom, its prob- 
ability of survival is calculated. 

v x 

Vy 

v z 

CM{\-ii(t,v r )) , 

GM(1-Mt,v r )) (3) 
GAf(l-/ifcv r )) . 

i ' 

r 2 

ufioir, t, 4>) 

In the equation above, /3o(r, t, <p) = r 2 yS(r, f, </>) is the scaled ion- 
ization rate at a heliographic latitude </>, r 2 — x 2 + y 2 + z 2 , GM is 
solar mass times constant of gravity, and v r {f) - (r(f) ■ v(t))/\r(t)\ 
is radial velocity of the atom. The position vector r(f) has coor- 
dinates x,y,z and velocity vector v(f) coordinates v x , v y , v z ; up is 
used to calculate the survival probability as 
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Vy 




v z 
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k U/3 ) 





exp(-w j8 ). 



(4) 



Once the position and velocity vectors of the test atom in the 
source region are found (respectively, r src , v src ), the value of the 
distribution function in the source region w src can be calculated. 
The value of the local distribution function is then be obtained as 
a product of w src and of the the survival probability Wj on . To ob- 
tain density and higher moments, the local distribution function 

is numerally integrated in the velocit y space d 3 V. 

The model of radiation pressure (Tarnopolski 2007) was de- 
veloped based on a series of observations of the solar Lyman- 
a line profiles, performed by SUMER/SOHO between solar 
minimum and maximum dLemaire et al.l |2002) and available 
on the Web. Since the SUMER instrument is located at the 
SOHO spacecraft orbiting a round the LI Lagr ange point be- 
tween the Earth and the Sun dWilhelm et al.| [l999). the observa- 
tions are free from geocoronal contamination that affected earlier 
measurements carri ed out from low-Earth orbits (e.g. OSO-5, 
IVidal-Madiatj|l975|). The absolu te calibration of the wavelength 
reported by Lemair e~et al.l d2002l) is better than +0.0015 nm. The 
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calibration of intensity was performed by direct comparison of 
results of integration of the observed profiles with the absolute 
fluxes obtained with the use of the SOLSTICE experiment. The 
accuracy of the net flux was reported at +10% level. It seems, 
however, that the data published on the Web are a smoothed 
version of actual measurements because despite the accuracy of 
10% reported in the paper the profiles do not show any scatter. 

Fittin g the func t ional forms of the profiles used in earlier 
studie s (lFahrlll978b IChabrillat & Kockartsl [19971: IScherer et alJ 
2000) did not yield sufficient accuracy of the model (cf Fig.Q]and 

To develop the model of n{v r , 7 t ot) discussed in this paper, the 
original data were rescaled to the spectral flux in the /v units as a 
function of radial velocity in km/s. A satisfactory re sult could be 
obtain ed with a 9-parameter function in the form ( Tarno polskil 
120071) : 

fi (v r , /tot(0) = A [1 + B 7 tot (f)] exp (-Cvfj (5) 

x [l + D exp (fv r - Gv 2 ) + H exp [-Pv r - gv 2 )] 

with the following parameters: 

A = 2.4543 x 10" 9 , B = 4.5694 x 10" 4 , C = 3.8312 x 10" 5 , 
D = 0.73879, F = 4.0396 x 10" 2 , G = 3.5135 x 10" 4 , 
// = 0.47817, P = 4.6841 x 10" 2 , Q = 3.3373 x 10" 4 . 
An interesting aspect of the model is that its sole solar- 
related parameter is the line-integrated flux in the Lyman-Q' 
line 7 to t. Any dependence of radiation pressure on time goes to 
the model via 7 to t. Hence, with a model of behavi o r of I tot in 
time on hand ( such a s, e.g., discussed by iBzowskH d200ll) and 
iBzowski et al.l d2008l) ) we can immediately construct a model of 
radiation pressure that will be dependent on radial velocity of 
the atoms as well. 

The function defined in Eq.((5]l with the parameters as listed 
above reproduce well the observed profiles both for solar min- 
imum (the profiles from July 27, 1996 to August 24, 1997) 
and for solar maximum (the profiles from August 20, 1999 to 
August 22, 2001). The accuracy of the fit for the most interest- 
ing region ±140 km/s about the line center exceeds the accuracy 
of observations declared by the authors. Furthermore, the line- 
integrated fluxes obtained from the model profiles agree well 
with the line-integra ted data and with the total fluxes reported by 
Tobiska et al. (2000) for a given day. The differences are on the 
order of 5%; only for the profile observed on May 20, 2000 the 
difference is 9%, which is still within the absolute calibration ac- 
curacy of 10%. A comparison of the experimental data with the 
fitted model is s hown in Fig.[T| A comparison of this model with 



fitted model is snown in rig.LU A comparison or tnis model witn 
earlier models dFahr et alJ ll98lb IChabrillat & Kockartsl Il997t 
IScherer et alJl2000l ) is presented in Fig. [2] 

As the absolute calibration of the line-integrated flux of the 
Sun changed since the time of publication, we show the origi- 
nal profiles recalibrated by multiplication by appropriate factor 
so that their absolute flux is in agreement with the values ac- 
cepted nowadays. The differences between the models in the ve- 
locity range +100 km/s around the line center are on the order of 
10 to 15%. For higher velocities the percentage differences be- 
tween the models are not so important because the absolute val- 
ues of radiation pressure are anyway small, and furthermore the 
fast atoms are n ot very sensitive to radiation pressure in general 
(Bzowski 2008). The highest differences in the absolute terms 
occur at the peaks of the profile, i.e., at the velocities charac- 
teristic for the fast wing of the distribution function of thermal 
populations of neutral hydrogen in the inner heliosphere. With as 
many as 9 data sets available, covering the whole span of solar 
activity level, we were able to come up with a model that seems 
to better reproduce physical reality than the former ones. 



This scheme was devised to enable calculating the density 
and higher moments of the distribution function in the case of 
a fully 3D and time dependent radiation pressure and ioniza- 
tion rate. The code has alread y been used in its full model by 
iTarnopolski & Bzowski (2008]) to assess the flux of inter stellar 
deuterium at the Earth orbit and bv lBzowski et al.l (|2008) to in- 
vestigate the density of neutral interstellar H at the termination 
shock based on the observations of hydrogen pickup ions by 
Ulysses. In the present study the code was restricted to the spher- 
ically symmetric and stationary case to facilitate assessment of 
the influence of the v, -dependence of radiation pressure on the 
distribution of interstellar hydrogen in the inner heliosphere, un- 
blended with other departures from the classical hot model. 



3. Calculations 

The simulations were performed in two groups. In the first one, 
we investigate the differences between the classical hot model 
and the Doppler model in a scenario with two populations at 
the termination shock: the primary interstellar population and 
the secondary population, which comes up up between the he- 
liopause and the bow shock from charge exchange between 
atoms from the original interstellar populations and the local 
heated and compressed plasma. These results are presented in 
Section 14.11 and will be referred to as the composite model, in 
the sense that it is composed of the two populations: primary 
and secondary. In the second group of simulations, we check the 
robustness of the conclusions drawn in Section 14.11 against un- 
certainties of the bulk velocity and temperature of interstellar gas 
at the termination shock and against uncertainties in the ioniza- 
tion rate in the inner heliosphere. The results of this series are 
discussed in Section |4~2l 

The calculations were performed on a mesh of heliocentric 
distances r equal to 0.4, 0.7, 1, 2, 4, 7, 10, 20, 50, and 100 AU 
and of offset angles from the upwind directions 6 equal to 0°, 
30°, 60°, 90°, 120°, 150°, and 180°, separately for solar mini- 
mum and solar maximum conditions, characterized by the fol- 
lowing values of the line-integrated solar Lyman-a' flux and net 
ionization rate: 

solar max.: / to t = 5.50 x 10 11 photons cirT 2 s _1 ,fl = 8 x 10" 7 s" 1 
solar min.: / tot = 3.53 x 10" photons crrT 2 s _1 ,/J = 5x 10" 7 s" 1 . 

The boundary conditions of the composite simulation were 
the parameters returned by the Moscow MC model at the nose 
of the termination shock for the following set of LIC conditions 
dlzmodenov et al.ll2003l) : 

proton density n Pt uc = 0.06 cm 4 , neutral H density «h,lic = 
0.18 cm" 3 , temperature of the gas Tuc = 6400 K, bulk veloc- 
ity v^lic = 26.4 km/s, the density of He + in the LIC « He + = 
0.008 cm" 3 . 

The resulting parameters of the primary and secondary pop- 
ulations at the nose of the termination shock were the following: 

Primary: «TS P ri = 0.03465 cm" 3 , v T s pr i = 28.512 km/s, 
r T s.p ri =6020K;' 

Secondary: «TSsec = 0.06021 cm" 3 , vxssec = 18.744 km/s, 
r TS ,sec = 16300 K. 

These parameters were adopted as boundary conditions for the 
composite simulations with the use of the Warsaw model. 

In addition to the dependence of radiation pressure on v r , the 
code inclu ded also a simplified m odel electron ionization rate, as 
adopted by Bzowsk i" et al.l (120081) . The model is spherically sym- 
metric and identical for solar minimum and solar maximum con- 
ditions. Inclusion of this model is a departure from the approach 
adopted in this study to keep all aspects of the modeling as close 
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Fig. 3. xg Classical/Doppler density ratios in the composite model for the solar minimum (left panel) and solar maximum conditions 
(right panel). The Sun is at the center in the (0,0) point, the upwind - downwind direction is along the y = line. The upwind, 
crosswind, and downwind axes are marked with heavy lines. 
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Fig. 4. xc Density excess of neutral interstellar hydrogen in the inner heliosphere, defined as a ratio of results of the classical hot 
model to the results of Doppler model. Left panel: density ratio at 1 AU as a function of offset angle from the upwind direction, 
right panel: density ratio at crosswind as a function of heliocentric distance. Solid lines correspond to solar minimum conditions, 
broken lines - to solar maximum conditions. 



to the classical hot model as possible. We decided to include the 
electron ionization to make sure that we do not miss possible 
disturbances of the Doppler-to-hot model ratios due to the extra 
ionization term which does not conform with the 1/r 2 distance 
profile of the charge exchange and photoionization rates. 

The simulations were repeated with the use of a model that 
differs from the former one only by the lack of sensitivity of radi- 
ation pressure to v r , i.e. effectively using the classical hot model 
with the electron ionization term added. The radiation pressure 
was equal to the spectral flux averaged by +30 km/s about 0, as 
indicated in FigQ] The numerical accuracy of the solutions was 
at the level of 1 to 2%. 

The simulations discussed in Section |4~2l were performed as- 
suming there is only one population of interstellar hydrogen at 
the termination shock. We explored the parameter space adopt- 
ing the net ionization rates at 1 AU /3, bulk velocities v src , and 
temperatures r src covering the entire range of expected values 
from the low values of (3 relevant for the fast solar wind con- 
ditions and solar minimum to the high values relevant for slow 
solar wind conditions and solar maximum, the ranges of v src and 
r src from the extremes predicted for the primary and secondary 
populations at the termination shock. 



This scan over f5 was performed assuming the gas temper- 
ature and velocity at the termination shock are equal, respec- 
tively, to 12500 K and 22 km/s, as derived bv lCosta et ail |l999) 
from an analysis of heliospheric Lyman-or glow in the upwind 
hemisphere, for 7 tot relevant for solar minimum and maximum 
conditions(see above). The scan over v src and r src space was per- 
formed assuming j3 and 7 tot values as listed for solar minimum 
and maximum conditions. 

The baseline set of simulations was performed with the use 
of the newly-developed Doppler model. Comparison simulations 
were performed with the use of the classical hot model on iden- 
tical spatial mesh and assuming identical parameters at the ter- 
mination shock and identical ionization rate. In this section the 
electron ionization was not included and hence the ionization 
rate was falling off as 1/r 2 . 

4. Results 

4. 1 . Composite model 

The simulations show that neglect of the v r -dependence of ra- 
diation pressure produces an excess of local hydrogen density 
(see Fig. 0). The excess q is defined here as the ratio of density 
returned by the classical hot model n c i ass to the density returned 
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Fig. 5. xd Density ratio of the classical hot model to Doppler 
model at the downwind axis as a function of heliocentric dis- 
tance. Solid line: solar minimum conditions, broken line - solar 
maximum conditions. 



by the Doppler model «dopp: q - «ciass/«dopp- Discrepancies be- 
tween the classical and Doppler models appear at ~ 3 AU and 
increase towards the Sun. Depending on the phase of the solar 
cycle and the offset angle from upwind, at 1 AU they are from a 
few dozen percent up to a factor of ~ 1 .7 (see Fig. [4). The excess 
density is a strong function of solar activity and while during so- 
lar minimum its typical values at 1 AU are about 25%, they are 
a little more than twice as large at solar maximum. Exception is 
a cone region around the downwind axis with the opening angle 
of ~ 60°, as shown in Fig. [4] and [5] at the downwind axis the 
hot model predicts a density deficit, which almost does not de- 
pend on the activity phase except at closest distances to the Sun. 
During solar maximum the excess is predicted even in the down- 
wind region, though only for distances smaller than ~ 1 .5 AU. 

The existence of the density excess can easily be understood 
after inspection of the shape of the solar Lyman-a- line profile, 
shown in Fig.Q] Since the profile has a minimum close to ra- 
dial velocity, the H atoms before approaching the Sun sense a 
stronger repulsion than in the case when one adopts a "flat" radi- 
ation pressure with the mean value averaged between +30 km/s 
around the line center. Hence in the Doppler case a larger por- 
tion of these atoms will be slowed down and repelled from the 
Sun, which reduces the density in comparison with the classi- 
cal case. The ensemble of surviving H atoms has a lower bulk 
velocity than predicted by the classical model, which leads to 
a further enhancement of the excess, because the ionization has 
more time to eliminate the atoms traveling with a smaller speed. 
On the other hand, there exists a sub-population of atoms which 
on approach towards the Sun experience a lower radiation pres- 
sure in the Doppler model than in the classical model. These are 
the atoms from the slow wing of the distribution function. In this 
case, although the Doppler model predicts less deceleration than 
the classical model, these atoms contribute comparatively little 
to the entire ensemble because they are more readily eliminated 
by ionization than the faster atoms. 

The Doppler model predicts different ratios of densities be- 
tween pairs of offset angles 6\, 62 (e.g. upwind/crosswind, up- 
wind/downwind etc.) than the predictions of the classical hot 
model. The discrepancies are not strong in the upper hemisphere, 
but escalate with the increase of difference 62 — 61. This must be 
one of the reasons for which interpretation of downwind obser- 
vations of the heliospheric Lyman-o- glow has always been more 
challenging than of th e upwind ones and usually returned dif- 
feren t conclusions (e.g. lLallement et aHll985at lOuemerai s et al.l 
1 1 992b . The challenge can be better appreciated when one realizes 
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Fig. 9. xk Classical/Doppler ratios of the fluxes of primary and 
secondary populations of neutral interstellar H atoms for solar 
minimum and maximum conditions as a function of offset angle 
from upwind. The lower pair of the lines corresponds to solar 
minimum, the upper pair to solar maximum. Within the solar 
minimum (lower) pair, the solid line corresponds to the primary 
population, the broken line to secondary population. Within the 
solar maximum (upper) pair, the dash-dot line is for primary, the 
dotted line for secondary population. 



that the wavelength-sensitivity effects are convolved with the ef- 
fects related to variations of solar Lyman-g flux and of the ion - 
ization rate, as discussed, e.g ., b vlRucinski & Bzowskil (fl995b); 
iBzowski & Ruciriskil (fl9933bT) : lBzowski et al.l (|2002|) . 

The Doppler and classical models predict surprisingly small 
differences between the local bulk velocities of the gas (Fig.|6]l. 
During solar minimum they are about 1 km/s and during solar 
maximum they increase to ~ 2 km/s. Differences in the radial 
component of the bulk velocities at 1 AU are also small (Fig. 
[7). They start only at the offset angle 6 60° and increase 
towards downwind; the classical model systematically predicts 
somewhat larger values. Since the differences in v r are small, we 
do not expect big differences in model spectral profiles of the 
heliospheric glow. 

Even though the differences in the local velocities are small, 
the effect of the wavelength dependence of radiation pressure on 
the magnitude of fluxes of the two heliospheric neutral H pop- 
ulations is not negligible because of the appreciable differences 
in the local densities of the gas. As shown in Fig.|9j during solar 
minimum the classical excess of the flux is 15 to 20%, depend- 
ing on the population, and during solar maximum it may reach a 
factor of 1 .4 to 1 .9 in the upwind hemisphere and even larger in 
the downwind region, but with the absolute magnitude of the flux 
reduced by 2 orders of magnitude from the solar minimum level 
(Fig. |8]l. In the upwind hemisphere it is very weakly sensitive to 
the offset angle. The wavelength dependence of radiation pres- 
sure must then be appropriately taken into account in interpreta- 
tion of direct in situ observations of neutral interstellar hydrogen 
atoms. The local flux of the atoms is a product of the local veloc- 
ity vector and of the local density, which is directly proportional 
to the density at the termination shock. Hence the density at the 
termination shock derived from an in situ-measured flux at 1 AU 
neglecting the effect of v r dependence of radiation pressure will 
be at least 15% off the mark (Fig|9]l, and most probably more, 
up to a factor of 2. 

Another factor that potentially can affect interpretation of 
such measurements are diff erences in the l ocal b ulk velocity vec- 
tors of the two populations. Mobi uTet al.l d2001) proposed to use 
the relative positions of the beams from the primary and sec- 
ondary populations to check whether their flow directions at the 
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Fig. 6. xl Bulk velocities of the primary and secondary populations of neutral interstellar H at 1 AU, calculated using the Doppler and 
classical models, shown as a function of offset angle from upwind for the solar minimum (left panel) and solar maximum conditions 
(right panel). Primary population: solid line (classical) and broken line (Doppler); secondary population: dotted line (Doppler) and 
dash-dot (classical). 
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Fig. 7. xi Radial velocities of the primary and secondary populations of neutral H atoms at 1 AU as a function of offset angle from 
upwind for solar minimum (left panel) and maximum conditions (right panel), calculated with the use of Doppler and classical hot 
model. The upper pairs of lines in both panels correspond to the primary population and the lower pairs to the secondary population. 
The classical hot model results are drawn, correspondingly, with solid and dash-dot lines, and the Doppler model results with broken 
and dotted lines. 
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Fig. 8. xj Absolute fluxes of neutral interstellar H atoms at 1 AU from the primary and secondary populations for solar minimum 
(left panel) and solar maximum conditions (right panel), shown as a function of the offset angle from upwind, calculated with the 
use of Doppler and classical hot models. The classical model results are drawn with solid line (primary populations) and dash-dotted 
line (secondary population); the Doppler model results with dashed (primary) and dotted lines (secondary). 



termination shock are parallel or not. Such measurements will 
hopefully b e soon realized by a for t hcom i ng NA SA SMEX mis- 
sion IBEX dMcComas et alj|2004 12005L 120061) . Fig. [TO] shows 
that the deflections of the two populations differ from the pre- 
dictions of the wavelength-independent model. For both popula- 
tions, in the region of the Earth orbit where detection by IBEX 
is the most probable the differences in deflection angle will be 



on the order of a few degrees, similar, but not identical for both 
populations. The magnitude of the deflection is a strong function 
of solar activity, with radiation pressure clearly playing the dom- 
inant role; the deflection itself is larger in the case of the Doppler 
model. 

Even though the differences are just a few degrees, they 
cannot be regarded as negligible. lOuemerais et al.l d 1999b and 
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Fig. 10. xm Deflection of the local upstream direction of the primary and secondary population of interstellar hydrogen at 1 AU as a 
function of offset angle from the gas inflow direction for solar minimum and maximum conditions. At solar minimum (left panel), 
less deflected is the primary population (solid line: classical model, broken line: Doppler model); the secondary is drawn with dotted 
line (Doppler) and dash-dotted line (classical). Same scheme for solar maximum plot (right panel). 



lLallement et al.l d2005l) suggested an upwind direction of neutral 
interstellar hydrogen that differs by about 5° from the inflow di- 
rection of interstellar helium measured bv lWittd d2004l) . The dif- 
ference is commonly attributed to a distortion of the heliospheric 
interface by the extraheliospheric magnetic field. Hence mostly 
affected should be the secondary population and the difference in 
the upwind direction of the primary and secondary populations 
inferred from analysis of in situ measurements interpreted with 
the use of a v r -insensitive model might be erroneously attributed 
to a deformation of the interface. Our results show that the de- 
flection of the secondary population obtained from the Doppler 
and classical models are very similar, but the deflections of the 
primary population differ by a few degrees. 

Finally, we assess the influence of the v,-sensitive radiation 
pressure on the expected PUI fluxes at 5 AU crosswind (i.e., 
the locatio n where the H + PUI fl ux observed by Ulysses is the 
strongest). Bzowski et al.l (120081) pointed out that the local pro- 
duction rate of pickup ions, which is a product of the local den- 
sity of neutral hydrogen and of the local ionization rate, is a 
weak function of radiation pressure and hence, necessarily, of 
its details including the v r dependence. They provided an esti- 
mate of the importance of the vv-dependence of radiation pres- 
sure for the local production rate of PUI which suggests it is 
on the level of a few percent. Here we report its influence of 
on the total flux of PUI, which was compute d fo llowing the sim- 
ple ap proach by Vasvliu nas & Siscoel (1 19761) and lGloeckler et al.l 
( 1993). Results of these simulations suggest that the classical hot 
model yields an excess flux with respect to the Doppler model, 
which at ~ 5 AU crosswind is equal to only ~ 7.5% at solar min- 
imum and to 8.5% at solar maximum. Thus we conclude that ne- 
glection of the v r dependence of radiation pressure results in an 
overestimate of the H + PUI flux at 5 AU from the Sun by about 
10%, irrespective of the solar activity phase. 

4.2. Comparison with the hot model: a scan of the parameter 
space 

The density excess q is a weakly increasing and almost linear 
function of the ionization rate (Fig.fTTTi. i.e. the higher the ioniza- 
tion rate one takes, the lower quality of the approximation given 
by the classical hot model one obtains. This can be explained 
by preferential eliminat ion of slower atoms from the local en- 
semble . As discussed bv lLallement et al.1dl985al) : lBzowski et al.l 
d!997l) . a result of ionization is a net acceleration of the gas by 



a few km/s because the atoms that have higher specific veloc- 
ities in the ensemble preferentially survive against ionization. 
Such atoms have also higher radial velocities and hence expe- 
rience a higher radiation pressure due to the non-flat shape of 
the solar Lyman-ff line profile. This force repels them from the 
Sun stronger than in the case when the solar line profile is flat. 
As a result, we have a reduction of density - consequently, the 
classical model predicts an excess of density, whose magnitude 
increases with the increase of the ionization rate. 

The sensitivity of the excess to the ionization rate increases 
with the offset angle from upwind and with the magnitude of 
the solar flux 7 tot . At the upwind axis for I tot corresponding to 
solar minimum conditions the slope dq/d/3 is equal to a half of 
the value relevant for crosswind and to ~ 0.1 of the value for 
downwind. For 7 tot relevant for solar maximum conditions the 
slopes are larger and the highest increase is for the upwind axis. 

In the downwind region, increasing the ionization rate im- 
proves the quality of the approximation provided by the classi- 
cal hot model with respect to the Doppler model. The parameter 
q is still an increasing function of the ionization rate, but since 
in the downwind region it is lower than 1, its increase means an 
improvement. 

The magnitude of the excess q is a weak function of temper- 
ature and bulk velocity of the atoms at the termination shock, as 
shown in Fig. [12] For 7 tot relevant for solar minimum conditions 
at 1 AU crosswind, it varies from 1.15 to 1.30 depending on the 
bulk velocity and temperature at the termination shock; for 7 to t 
relevant for solar maximum conditions the amplitude is similar, 
but the value of q larger. 

Generally, the quality of the hot model approximation of 
density is consistently lower for higher values of / to t, which is 
modulated by solar activity level. The sensitivity of radiation 
pressure to radial velocity thus modifies the relations between 
densities during various phases of solar activity (the amplitude o f 
modulation), discussed, e.g. , by iRucinski & Bzowski (1995b); 
iBzowski & Ruciriskil d!995al) and lBzowski et al l (120021) . 

On the other hand, it seems that the effect on the model 
values of the heliospheric Lyman-a glow should be on the or- 
der of 10%. An exact calculation requires developing a model, 
where the local bulk velocity and temperature of the gas and 
their radial gradients will be taken into account, together with 
the wavelength dependence of the solar illuminating flux and its 
wavelength-differential attenuation, which increases with the in- 
crease of the solar distance. We believe this effect to be relatively 
small because during solar minimum, when the contribution to 
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Fig. 11. xf Classical hot/Doppler model density ratios for the upwind (dots), crosswind (dashed) and downwind axes (solid) at 1 AU 
as a function of the ionization rate for solar minimum (left panel) and maximum conditions (right panel). 
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Fig. 12. xe Density excess of the classical hot model with respect to the Doppler model of neutral interstellar hydrogen as a function 
of temperature and bulk speed in the source region for solar minimum (left panel) and solar maximum conditions (right panel). 



the net backscatter intensity from the gas near the Sun is largest, 
the q factor is lowest. During solar maximum, when q is largest 
and hence the quality of the v r -independent model is lowest, the 
local density of the gas is very much reduced, so increased is the 
contribution to the net signal of more distant regions of the helio- 
sphere, where the q values are relatively small. In consequence, 
the v,. dependence of radiation pressure should appreciably affect 
in situ measurements within ~ 5 AU from the Sun, but only in 
a small degree the photometric observations of the heliospheric 
glow. 

5. Conclusions 

We performed a comparison of density, bulk velocity and flux of 
neutral interstellar hydrogen in the inner heliosphere, calculated 
using either the classical hot model or a newly-developed model 
with radiation pressure being a function of specific radial veloc- 
ities of individual atoms. The conclusions are the following. 

1. Differences between the Doppler and classical hot models 
are restricted to ~ 5 AU from the Sun. The classical hot 
model overestimates the density of neutral hydrogen gas 
everywhere in the inner heliosphere except a ~ 60° cone 
around the downwind direction, where a density deficit is 
predicted. Generally, the density excess/deficit is a strong 
function of the offset angle from the upwind direction and 
of the heliocentric distance (Figs[3]|4]and[5]l. 



2. The magnitude of density excess varies with the solar ac- 
tivity level and is higher at solar maximum (Figs [3] and 
@J, i.e. v, -independent models overestimate the local density 
stronger during solar maximum. 

3. The density excess is a weakly increasing, almost linear 
function of the ionization rate in the inner heliosphere and a 
weak almost bi-linear function of the bulk velocity and tem- 
perature of the gas at the termination shock (Figs[TTIand[T2l. 

4. The classical model shows a higher upwind-downwind am- 
plitude of radial velocities at 1 AU than the Doppler model; 
the differences are on the order of 10%. The excess of ab- 
solute flux of neutral atoms returned by the classical model 
are on the order of 15% at 1 AU during solar minimum, i.e. 
at the time of observations by the forthcoming IBEX mis- 
sion. The excess is much higher during solar minimum, but 
the absolute values of the fluxes are lower by a few orders of 
magnitude than during solar minimum (Figs |7|8l9l6t . 

5. The classical model returns a different deflection of the local 
upstream directions of the primary population of interstel- 
lar neutral H than the Doppler model. The deflections of the 
secondary population returned by the classical and Doppler 
model are practically identical. Hence interpretation of in 
situ measurements of these fluxes must be performed very 
carefully, with appropriate account for the solar line profile 
shape to avoid confusing deflections caused by the Doppler 
effect in radiation pressure with, e.g., results of deformation 
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of the heliospheric interface due to external magnetic field 

(Fig. Hi). 

6. Discrepancies between the PUI flux at Ulysses in aphelion 
returned by the Doppler and classical models are on the order 
of 10% and are almost insensitive to the phase of solar cycle. 
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